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This article is a pedagogical review of Monte Carlo methods for the self-avoiding walk, with emphasis on the 
extraordinarily efficient algorithms developed over the past decade. Many more details can be found in [1]. 



1. INTRODUCTION 

This talk has no direct relation to QCD; it's 
therefore intended as entertainment. My goal is 
to show a nontrivial statistical-mechanics prob- 
lem for which we have successfully developed 
collective-mode algorithms that completely elimi- 
nate the critical slowing-down. Hopefully some of 
these ideas can be adapted to problems in random 
surfaces, and possibly even to spin and gauge the- 
ories. Anyone who wants more details (and more 
references) is referred to my review article [1]. 

So, what is a self-avoiding walk (SAW)? Obvi- 
ously, it's a nearest-neighbor path (on some given 
lattice) that never visits any site more than once. 
Figure 1 shows some typical TV-step SAWs with 
N = 100, 1000, 10000 on the square lattice, scaled 
down by a factor N'^ with i/ = | (that turns 
out to be the right critical exponent in dimension 
d = 2). Already we can see that a long SAW is 
going to be a fractal object. Moreover, the shape 
and size of these SAWs are completely different 
from the shape and size of ordinary random walks 
(which have i/ = |); this is a different universality 
class. Note, finally, that even by looking closely 
at the 10000-step SAW one can't tell whether it's 
strictly self-avoiding or merely self-repelling; nev- 
ertheless, the universality class doesn't depend on 
such short-distance details. 

The purpose of this talk is to explain to you 
how I generated these three SAWs. 

But first of all one should ask: Who cares? 
Why is the SAW an interesting model? The 
answer is twofold. On the one hand, the SAW 
is a model of a linear polymer molecule (such 
as polystyrene) in a good solvent. Further- 
more, long chains are experimentally accessible: 




Figure 1. Typical TV-step SAWs with N = 
100, 1000, 10000 on the square lattice, scaled 
down by a factor N" with = 3/4. The 10000- 
step SAW was not prepared specially for this con- 
ference in Australia. 
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one can reach N as high as ~ 10^. Secondly, 
the SAW is isomorphic to the 0(n) cr-model [also 
known as the n-vector model] analytically contin- 
ued to n = 0; more precisely, the SAWs are the 
strong-coupling diagrams of the cr-model. So the 
SAW belongs to a family we already know and 
love: its siblings include the Ising (n = 1), XY 
(n = 2), classical Heisenberg (n = 3) and spheri- 
cal (n = oo) models. 

The long-chain limit TV ^ oo of the SAW is 
thus a critical phenomenon in the usual sense: 
it's described by a continuum quantum field the- 
ory, so that one expects universal behavior with 
the usual plethora of critical exponents, universal 
amplitude ratios, universal scaling functions, and 
so forth. Many of these theoretical predictions 
are experimentally testable by light scattering on 
dilute polymer solutions. 

In fact, the SAW is an exceptionally favorable 
"laboratory" for the numerical study of critical 
phenomena, for three reasons: 

1) There are no finite-volume effects; one can 
study directly an TV-step SAW in an infinite lat- 
tice. 

2) There is no L"^ factor in the computational 
work; one doesn't have to worry about the sites 
where the walk might be but isn't. 

3) There is, to be sure, critical slowing-down; 
but this can be overcome by cleverly designed 
collective-mode algorithms. 

2. DEFINITIONS 

Let cjv [resp. cjv(a;)] be the number of TV-step 
SAWs on Z** starting at the origin and ending 
anywhere [resp. ending at x]. These quantities 
are believed to have the asymptotic behavior 

CM ~ ti'^N^-' (1) 
CNix) ~ fj,^ N""-'"-^ {x fixed ^ 0) (2) 

as TV ^ oo. [Compare ordinary random walks, 
for which CAT = (2d)^ and c^ix) ~ {2d)^N-'^'\] 
Here ji is called the connective constant of the lat- 
tice: it is like a critical temperature, and there- 
fore lattice-dependent. By contrast, 7 and aging 
are critical exponents, and thus expected to be 
universal among lattices of a given dimension d. 
Now consider, for each fixed TV, the probability 



distribution in which each TV-step SAW gets equal 
weight. Then we can ask about the mean-square 
end-to-end distance (R^), the mean-square radius 
of gyration (Rg), and so forth: all these measures 
of the size of the SAW are believed to grow as 

{RI)n, {RI)n ~ TV^-^ (3) 

as TV ^ 00, where v is another (universal) criti- 
cal exponent. Moreover, the ratio (R^) / (R^) — 
along with many other dimensionless measures of 
the shape of a SAW — tends as TV ^ 00 to a 
universal limiting value. 

More information on the SAW can be found in 
the excellent book by Madras and Slade [2]. 

Remark: In dimension d = 2, some of these 
universal quantities have been computed exactly 
by conformal-field-theory and Coulomb-gas tech- 
niques: for example, v = ^ and 7 = || [3]. It's 
an open question whether all of them can be so 
computed; for example, the exact value of the lim- 
iting ratio (i2^)/(i2g) is unknown (to four decimal 
places it's 0.1403 [4-6]). 

Different aspects of the SAW can be probed in 
four different ensembles: 

• Fixed-length, fixed-endpoint ensemble 
(fixed TV, fixed x) 

• Fixed-length, free-endpoint ensemble (fixed 
TV, variable x) 

• Variable-length, fixed-endpoint ensemble 
(variable TV, fixed x) 

• Variable-length, free-endpoint ensemble 
(variable TV, variable x) 

The fixed-length ensembles are best suited for 
studying the critical exponent v, while the 
variable-length ensembles are best suited for 
studying the connective constant ji and the crit- 
ical exponents aging (fixed-endpoint) or 7 (free- 
endpoint). Physically, the free-endpoint ensem- 
bles correspond to linear polymers, while the 
fixed-endpoint ensembles with |a;| = 1 correspond 
to ring polymers. All these ensembles give equal 
weight to all walks of a given length; but the 
variable-length ensembles have considerable free- 
dom in choosing the relative weights of different 
chain lengths TV. 
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Most of the algorithms I'll discuss here work 
in the fixed-TV, free-endpoint ensemble. Indeed, 
somewhat less progress has been made in elimi- 
nating critical slowing-down in the other three en- 
sembles, and many fascinating problems remain 
open. See [1] for details. 

3. STATIC ALGORITHMS 

I call a Monte Carlo algorithm static if it pro- 
duces a sequence of statistically indepdendent 
samples: that is, it's of the form "call the subrou- 
tine and it returns a random SAW" , independent 
of all previous ones. 

The most obvious static technique for gener- 
ating a random TV-step SAW is simple sampling: 
just generate a random TV-step ordinary random 
walk, and reject it if it is not self-avoiding; keep 
trying until success. It is easy to see that this 
algorithm produces each TV-step SAW with equal 
probability. Of course, to save time we should 
check the self-avoidance as we go along, and reject 
the walk as soon as a self-intersection is detected. 

The trouble with this algorithm is, of course, 
the exponentially rapid sample attrition for long 
walks: the probability of an TV-step walk be- 
ing self-avoiding is CN/[2d)^ ~ [fi/2d)^. Some 
improvement can be obtained by modifying the 
walk-generation process so as to produce only 
walks without immediate reversals; but the suc- 
cess probability still decays like [fi/[2d — 1))^ ■ 
One can try higher-order variants of simple sam- 
pling, in which walks with loops of length < r 
are automatically absent; but now it's non-trivial 
to ensure that each TV-step SAW gets generated 
with equal probability, and one still has exponen- 
tial attrition (albeit a weaker one). All in all, this 
approach seems to be a dead end. 

So it's perhaps surprising that there is a rea- 
sonably efficient static algorithm for generating 
SAWs: it's called dimerization, and it's an im- 
plementation of the computer scientists' princi- 
ple of "divide and conquer". To generate an TV- 
step SAW, we generate two independent (TV/2)- 
step SAWs ("dimers") and attempt to concate- 
nate them. If the result is self-avoiding, we are 
done; otherwise, we discard the two walks and 
start again from scratch. This procedure can now 



be repeated recursively: to generate each of the 
(TV/2)-step SAWs, we generate a pair of (TV/4)- 
step SAWs and attempt to concatenate them, and 
so on. For TV < some cutoff TVq, we stop the re- 
cursion and generate the SAWs by some prim- 
itive method, such as non-reversal simple sam- 
pling. The dimerization algorithm can thus be 
written recursively as follows: 

function dim(TV) 

if TV < TVo then 

u! <— nrssamp(TV) 

return w 
else 

TVi ^ [TV/2J (integer part) 
TV2 ^ TV - TVi 
start: 
^^(1) ^ dim(TVi) 



dim(TV2) 



w <— a'(^) o ai(2) (concatenation) 
if w is not self-avoiding goto start 
return w 
endif 

It is easy to prove inductively that algorithm dim 
produces each TV-step SAW with equal probabil- 
ity, using the fact that the subroutine nrssamp 
does so. It is crucial here that after a failure we 
discard both walks and start again from scratch. 

Let's analyze the efficiency of the dimerization 
algorithm under the scaling hypothesis 



CAT ^ Afi^N^-'^ . 



(4) 



Let Tjv be the mean CPU time needed to gener- 
ate an TV-step SAW by algorithm dim. Now, the 
probability that the concatenation of two random 
(TV/2)-step SAWs yields an TV-step SAW is 



Pn 



(cJv/2)^ 



(5) 



where B = A/4:'^~^. We will need to generate, on 
average, 1/pjv pairs of (TV/2)-step SAWs in order 
to get a single TV-step SAW; hence 



l-N 



BTVT-^2T, 



N/2 



(6) 



(We have neglected here the time needed for 
checking the intersections of the two dimers; this 
time is linear in TV, which, as will be seen shortly. 
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is negligible compared to the time 2Tjv/2 for gen- 
erating the two dimers.) Iterating this k times, 
where k = log2{N / Nq) is the number of levels, we 
obtain 

where 

Ci = (t-1)/2 (8) 
C2 = (5 - 37)/2 + log2 A (9) 

and Co depends on Nq- Thus, the growth of Tjv is 
slower than exponential in N; but if 7 > 1 (which 
occurs for d < 4) it is faster than any polynomial 
in N . Fortunately, however, the constants Ci and 
C2 are very small, so that in practice Tjv behaves 
like N~^~^ up to N of order several thousand 
(resp. several million) in d = 2 (resp. d = 3). 

This may be the only known subexponential- 
time static algorithm for a nontrivial statistical- 
mechanical problem.^ It's an open question 
whether there exists a polynomial-time algorithm 
for generating SAWs in dimension d < 4. 

4. DYNAMIC ALGORITHMS 

Most of the Monte Carlo algorithms familiar to 
quantum field theorists are dynamic algorithms: 
that is, they generate a sequence of correlated 
samples from some Markov process having the de- 
sired probability distribution as its unique equi- 
librium distribution. The main problem with 
dynamic algorithms, as we all know, is critical 
slowing- down^: the autocorrelation time r of the 
Monte Carlo dynamics^ grows as the critical point 
is approached. For the SAW, "criticality" means 
TV" ^ 00. 

The elementary moves in a SAW Monte Carlo 
algorithm can be classified according to whether 
they are 

^ Here "nontrivial" is meant to exclude models like inde- 
pendent percolation, Gaussian fields, etc. 
^ See [7] for an introduction to dynamic Monte Carlo meth- 
ods and critical slowing-down. 

^In fact there are several distinct autocorrelation times — 
notably the exponential autocorrelation time Texp and the 
integrated autocorrelation times Tj^^ ^ for various observ- 
ables A — and (contrary to much belief) these can have 
different dynamic critical exponents Zexp and Zi^t ^ . See 
[8] for discussion. 



• TV-conserving or TV-changing 

• endpoint-conserving or endpoint-changing 

• local, bilocal or non-local 

Obviously, fixed- TV algorithms must use only TV- 
conserving moves, while variable-TV algorithms 
are free to use both TV-conserving and TV- 
changing moves (and indeed must use some of the 
latter in order to satisfy ergodicity). An anal- 
ogous statement holds for fixed- a; and variable- 
X algorithms with regard to endpoint-conserving 
and endpoint-changing moves. 

The most important distinction is between lo- 
cal, bilocal and non-local moves. Pure local 
and bilocal algorithms are easy to invent, but 
they lead to critical slowing-down (just as for the 
analogous algorithms in QCD). Non-local algo- 
rithms are harder to invent, but they offer at least 
the possibility of radically reducing or even com- 
pletely eliminating the critical slowing-down. 

Obviously there's no space here to review in 
detail all the known algorithms for the SAW (see 
[1]). Instead, I'd like to present one algorithm 
of each type, just to give the flavor of how these 
algorithms work and how they can be analyzed. 

4.1. Local algorithms 

A local move is one that alters only a few con- 
secutive sites ("beads") of the SAW, leaving the 
other sites unchanged. Otherwise put, a local 
move excises a small piece from the original SAW 
and splices in a new local configuration in its 
place. (Of course, it is always necessary to ver- 
ify that the proposed new walk is indeed self- 
avoiding.) Let's restrict attention for simplicity 
to N -conserving local moves. 

Figure 2 shows all the possible one-bead local 
moves (on a hypercubic lattice). Move A is a 
"one-bead flip" (also called "kink-jump"); it is 
the only one-bead internal move. Moves B and C 
are end-bond rotations. 

Figure 3 shows all the possible internal two- 
bead moves. Move D is a "180° crankshaft". 
Move E is a "90° crankshaft"; of course it is pos- 
sible only in dimension d > 3. Move F is a "two- 
bead L-flip". Move G permutes three successive 
mutually perpendicular steps (which lie along the 
edges of a cube); again this is possible only in di- 
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Figure 2. All one-bead local moves. (A) One-bead flip. (B) 90° end-bond rotation. (C) 180° end-bond 
rotation. 
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Figure 3. All internal two-bead local moves. (D) 180° crankshaft. (E) 90° crankshaft (d > 3 only). (F) 
Two-bead L-flip. (G) Cube permutation (d > 3 only). 



mension d > 3. I leave it to the reader to con- 
struct the list of two-bead end-group moves. 

A local algorithm can be built by combining 
any subset of local moves: for example, one pop- 
ular algorithm (used heavily by chemical physi- 
cists) employs moves A-E. Unfortunately, all lo- 
cal TV-conserving algorithms have a fatal flaw: 
they are nonergodic for sufficiently large N.'^ For 
algorithms based on moves of k or fewer beads, 
nonergodicity arises in dimension d = 2 for all 
N > 16k + 63, and for quite a few smaller N as 
well [9, Theorem 1]. 

Even if we put aside the problem of noner- 
godicity, the local algorithms have severe crit- 
ical slowing-down. A plausible heuristic argu- 
ment suggests that r ~ 7\r2+2i/ _ ^>3. 
sider the evolution of the center-of-mass vector 
of the chain. At each elementary move, a few 
beads of the chain move a distance of order 1; so 
ArcM ~ 1/-/V. But in order to traverse its equilib- 



rium distribution, fcm must change by something 
of order its standard deviation, which is ~ N'^ . 
Assuming that fcm executes a random walk, this 
takes a time N^'^^'^ . It would be interesting 
to know whether this prediction for the dynamic 
critical exponent is exact or merely approximate.^ 
In any case, the N-^ critical slowing-down of the 
local algorithms makes it difficult in practice to 
get beyond N ~ 100. 

4.2. Bilocal algorithms 

A bilocal move is one that alters two disjoint 
small groups of consecutive sites (or steps) of the 
walk; these two groups may in general be very far 
from each other. Here are some examples: 

• The slithering- snake (or reptation) move, 
which deletes a bond from one end of the walk 
and appends a new bond (in an arbitrary direc- 
tion) at the other end [Figure 4]. 

• The kink transport move, which deletes a kink 



^At least in dimensions d = 2,3. The result is probably 
true also in dimensions d > 4, but it has not yet been 
proven. 



^Usually dynamic critical exponents are not expressible 
in terms of static ones, except for trivial (Gaussian-like) 
models. 
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Figure 4. The slithering-snake (reptation) move. 
Head of the walk is indicated by x . Dashed lines 
indicate the proposed new step (resp. the just- 
abandoned old step). 
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Figure 5. The kink-transport move. A kink has 
been cleaved from AB and attached at CD. Note 
that the new kink is permitted to occupy one or 
both of the sites abandoned by the old kink. 
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Figure 6. The kink-end reptation ( >) and end- 
kink reptation (< — ) moves. In >, a kink has 

been cleaved from AB and two new steps have 
been attached at the end marked x . Note that 
the new end steps are permitted to occupy one or 
both of the sites abandoned by the kink. 



at one location along the walk and inserts a kink 
(in an arbitrary orientation) at another location 
[Figure 5]. 

• The kink-end reptation move, which deletes a 
kink at one location along the walk and appends 
two new bonds (in arbitrary directions) at one of 
the ends of the walk [Figure 6 — >]. 

• The end-kink reptation move, which deletes 
two bonds from one of the ends of the walk and in- 
serts a kink (in an arbitrary orientation) at some 
location along the walk [Figure 6< — ]. 



Bilocal algorithms can be either nonergodic 
(e.g. pure reptation®) or ergodic (e.g. various 
combinations). As for the critical slowing-down, 
a plausible heuristic argument suggests that r ~ 
N^: the SAW transforms itself by random back- 
and-forth slithering along the chain; after 
moves, this slithering will have carried it N steps, 
at which time all the original bonds of the chain 
will have disappeared and been replaced by ran- 
dom new ones. It would be interesting to know 
whether this prediction for the dynamic critical 
exponent is exact or merely approximate. In 
practice, with a algorithm one can reach 

N ~ 1000. 

The study of bilocal algorithms is still in its 
infancy. We need to understand better the issues 
of ergodicity and critical slowing-down, and the 
delineation of dynamic universality classes. 

Remark: The BS algorithm [10] is a variable- 
N algorithm closely related to reptation: it is er- 
godic and can be proven to have r ~ N~^. 

4.3. Non-local algorithms 

The possibilities for non-local moves are almost 
endless, but it is very difficult to find one which 
is useful in a Monte Carlo algorithm. There are 
two reasons for this: Firstly, since a non-local 
move is very radical, the proposed new walk usu- 
ally violates the self-avoidance constraint. (If you 
move a large number of beads around, it becomes 
very likely that somewhere along the walk two 
beads will collide.) It is therefore a non-trivial 
problem to invent a non-local move whose accep- 
tance probability does not go to zero too rapidly 
as TV ^ oo. Secondly, a non-local move usually 
costs a CPU time of order N (or in any case 
with p > 0), in contrast to order 1 for a local or 
bilocal move. It is non-trivial to find moves whose 
effects justify this expenditure (by reducing Tint, a 
more than they increase Tcpu)- 

One extremely successful non-local algorithm 
is the pivot algorithm [11-13]. Here we choose 
some site along the walk as a pivot point, and ap- 
ply some symmetry operation of the lattice (e.g. 
rotation or reflection) to the part of the walk sub- 
sequent to the pivot point [Figure 7]. It's easy to 



Exercise for the reader: Find a frozen configuration. One 
solution can be found in [1, Section 6.4.2]. 
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Figure 7. The pivot move (here a +90° rotation). 
The pivot site is indicated with an x . Dashed 
lines indicate the proposed new segment. 



prove detailed balance, and with some work er- 
godicity can be proven as well [13]. 

At first thought this seems to be a terrible 
algorithm'^: for N large, nearly all the proposed 
moves will get rejected. In fact, this latter state- 
ment is true, but the hasty conclusion drawn from 
it is radically false! The acceptance fraction / 
does indeed go to zero as TV ^ oo, roughly like 
N~P; empirically, it is found that the exponent 
p is ^ 0.19 in d = 2 and ^ 0.11 in d = 3. But 
this means that roughly once every moves 
one gets an acceptance. And the pivot moves are 
very radical: one might surmise that after very 
few accepted moves (say, 5 or 10) the SAW will 
have reached an "essentially new" configuration. 
One conjectures, therefore, that the autocorre- 
lation time r of the pivot algorithm behaves as 
On the other hand, a careful analysis of 
the computational complexity of the pivot algo- 
rithm shows that one accepted move can be pro- 
duced in a computer time of order N.^ Combin- 
ing these two facts, we conclude that one "effec- 
tively independent" sample (at least as regards 
global observables) can be produced in a com- 
puter time of order N. This is vastly better than 



^Indeed, this was my initial view. See [13, footnote 3] for 
a mea culpa. 

^Things are in fact somewhat more subtle: heuristic argu- 
ments and numerical evidence suggest that Tj^^ ^ ~ 
for observables A that are "essentially global" (like i?g and 
i?^), while Tint, A ^ N^'^^ for "local" observables (like the 
number of 90° angles in the walk). See [13, Sections 3.2, 
3.3, 4.2 and 4.3] for details. 

^One should use a hash table to test self-avoidance, and 
look for a self-intersection starting at the pivot point and 
working outwards. Failed moves can then be detected 
in a mean CPU time of order N^~^ (because the self- 
intersection is most likely to be near the pivot point), while 
successes of course require a time of order N . See [13, Sec- 
tions 3.4 and 4.4] for details. 
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Figure 8. Interpenetration ratio ^ versus chain 
length N, for SAWs in d = 3. Error bar is one 
standard deviation. Data from [5]. 



the 7\r2+2i/ q£ ^Yie local TV-conserving algorithms 
and the of the bilocal algorithms. Indeed, 

this order of efficiency cannot be surpassed by any 
algorithm which computes each site on successive 
SAWs, for it takes a time of order N simply to 
write down an TV-step walk! 

In practice, with the pivot algorithm one can 
reach chain lengths N of order 10^, with high 
statistics [5]. One can then obtain extremely 
high-precision (±0.001) estimates of critical expo- 
nents and universal amplitude ratios, with good 
control over corrections to scaling}'^ An example 
is shown in Figure 8, which plots the interpene- 
tration ratio ^ for pairs of TV-step SAWs (it's a 
kind of dimensionless renormalized coupling con- 
stant) as a function of TV. We can estimate the 
limiting value = 0.2471 ± 0.0003, and we can 
see clearly the strong corrections to scaling. This 
plot also shows that most of the last 40 years 
of polymer theory is wrong (albeit fixable), but 
that's another story [14,5]. 



^^Corrections to scaling are a very serious problem when 
one is aiming for this level of accuracy. For example, in the 
d = 3 SAW, to get the exponent u with a systematic error 
of less than 0.001 from a pure power-law fit (J?^) ~ N^^ , 
one must throw away all the walks with iV < 1000 [5]! 
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4.4. Other algorithms 

Well, there are many more interesting algo- 
rithms that I don't have time to discuss. For 
the variable- TV, free-endpoint ensemble we have 
the slithering-tortoise (BS) [10] and the join-and- 
cut [15] algorithms. For the variable- TV, fixed- 
endpoint ensemble we have the BFACF [16] and 
BFACF + cut-and-paste [17] algorithms. But 
there's lots of room for inventing new and better 
non-local algorithms for these ensembles. And for 
all the algorithms described here, we need to un- 
derstand better their dynamic critical behavior. 

A big open problem is to devise good non-local 
algorithms for the SAW with nearest-neighbor 
self- attraction, which has a collapse transition 
(called the theta point) as the temperature is var- 
ied. The pivot algorithm (for example) can eas- 
ily be modified to handle a nearest-neighbor in- 
teraction, by inserting a Metropolis accept/reject 
step; but its efficiency deteriorates markedly in 
the neighborhood of the theta point and even 
more drastically in the collapsed regime. 

Another open problem is to generalize these al- 
gorithms to related models such as branched poly- 
mers, random surfaces, etc. For these problems 
local algorithms are known, but good non-local 
algorithms are just beginning to be invented. 

5. CONCLUSIONS 

The SAW is an extraordinary success story for 
the development of non-local (collective-mode) al- 
gorithms: the CPU time needed to generate one 
"effectively independent" SAW has been reduced 
from order e^-^ (simple sampling) to TV"^' i°g^+C2 
(dimerization) to N-^ (local algorithms) to 
(bilocal algorithms) to N (pivot algorithm). 

Here's the recipe for success of a non-local algo- 
rithm: Offer the system a collective-mode move 
that it wants. (The move has to be radical but 
sensible.) And do it in a CPU time of order N 
(~ volume in a spin or gauge system) or in any 
case not much more. Of course, this is not really 
a "recipe"; it's more like "guidelines". To invent 
a good non-local algorithm, one needs physical 
insight combined with cleverness. 

The big open question is: Can any of these al- 
gorithms from polymer physics inspire successful 



collective-mode algorithms for spin or gauge the- 
ories? Or vice versa? 

Oh, yes, I almost forgot: the three SAWs in 
Figure 1 were produced using the pivot algorithm. 

This research was supported by the Petroleum 
Research Fund and by NSF grant DM S- 9200719. 
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